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ABSTRACT 

We consider the decelerating shock instabihty of Vishniac for a finite layer of con- 
stant density. This serves both to clarify which aspects of the Vishniac instability mech- 
anism depend on compressible effects away from the shock front and also to incorporate 
additional effects of finite layer thickness. This work has implications for experiments 
attempting to reproduce the essential physics of astrophysical shocks, in particular their 
minimum necessary lateral dimensions to contain all the relevant dynamics. 

Subject headings: hydrodynamics - instabilities - shock waves 



1. Introduction 



Vishniac ( 1983 ) outlined a theory of instabilities for a system of a decelerating shock accreting 



mass, modeled as a thin mass shell layer possessing no internal structure. In [Vishniac and Ryu 



( 1989 ), the theory was expanded to include a layer of post-shock material, exponentially attenuating 



in density. Other work ( Bertschinger (1986); Kushnir et al. (2005); among others) has described 
the perturbation of self-similar solutions for the post-shock flow. The present work complements 
these investigations, modeling the post-shock flow as a finite thickness layer of constant density and 
considering both compressible and incompressible post-shock states. This allows us both to more 
clearly understand which mechanisms depend on the compressibility of the shocked gas and which 
are common to any shock system undergoing deceleration. 

Early in the lifetime of an impulsively driven shock, when the post-shock layer thickness is 
small compared to its compressible length scale, an exponential scale cannot be formed and the 
density profile may be closely approximated by a square wave, as a fluid everywhere of constant 
density. In the shock's frame, upstream fluid is entering the shock with a speed Vg and exiting 
it with a speed U = Vgrj, where ry is the inverse compression ratio associated with the shock, 
including both the viscous density increase and any subsequent, localized further density increase 



in consequence of radiative cooling ( Drake 2006 ) . 
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Fig. 1. — Schematic of the decelerating shock system. The sohd black line is the shock, the dashed 
line above the dense rear layer is the rear material interface. The left-hand arrow depicts the 
in- frame inertial force with acceleration (— V^). 



We have also in mind throughout this paper experiments (Reighard et al. 2006: Bouquet 



et al.|[2004 Bozier et al.||1986 ) that have been carried out to investigate radiating shock dynamics. 
Experiments of this type are often designed to be scaled to relevant astrophysical investigation 
( Remington et al.|2006 ) . The particular experiments of Reighard et al. ( 2006 ) feature characteristic 
shock velocities of over 100 km/sec, shock tubes with 625 /im diameters, and strong deceleration 
throughout the shocks' lifetimes. 



2. System of a Decelerating Shock with a Dense Downstream Layer 

We consider the shock in its own, decelerating frame. The system is depicted in Figure [l] The 
shock is placed at z = 0, with flow entering it from the negative z direction at speed Vg, density 
pO) and with negligible thermal pressure. Flow is exiting the shock toward positive z with speed 
U , density p, and isotropic pressure P. We will model the downstream, rear layer as a constant 
density region of finite, increasing thickness from z = Q lo z = H. The rear surface of the dense 
layer will be taken to be a free interface at constant pressure. Beyond the rear layer will be taken 
as a region of constant thermal pressure Pi. 

The native surface wave modes in the system will be right- and left-propagating waves on the 
two surfaces of the dense layer, leading to four modes in total. As drawn in Figure [T| the upper 
surface, a material discontinuity, is stable if the shock frame is decelerating and is characterized 
by surface gravity modes. The bottom surface, a shock at which compressibility is not suppressed, 
will feature propagating acoustic modes. The waves which appear in our coupled system will be 
modifications of these waves which appear on these surfaces in isolation. In particular, the modified 
acoustic waves along the shock surface will be identified as bending modes of the entire dense layer. 

In order to understand the fundamental cause of the instability, we will here be considering 
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the fluid both ahead and behind the shock to be held at (different) densities constant in both 



space and time. This practice is described and defended by Hayes and Probstein (1966), who 
in their book "consider constant-density hypersonic flows, though we should never consider the 
fluid in a hypersonic flow as incompressible." The pressure profile behind the shock is hydrostatic, 
P{z) = Pi — {H — z)pVs, which leads to increasing pressure at the shock front when the shock is 
decelerating. Perturbations to density by the waves under investigation will be discussed. 



3. Linear Perturbations of the System 
3.1. Solutions Inside the Post-Shock Fluid 

We begin with the inviscid fluid equations 

p{dtv + v-Vv) = -VP - pVsZ (1) 
dtp + v-Vp = -pV • V (2) 

with total velocity v = {u,0,w + U) and P = P + 5P. We will insert the perturbation 6p only 
in the continuity equation; the coupling of 6p to the frame's acceleration will be suppressed. This 
allows us to ignore mode purely internal to the layer, concentrating on the overall shock and layer 
system. Since log p/po » log {p + Sp)/p for any reasonable density perturbations, we expect the 
dynamics of the system to be dominated by the compression at the shock. The omission of the 
term dpVg is also required for consistency with the assumption of our square wave density profile; 
the system will otherwise begin to evolve into an exponential atmosphere. 

We first let the perturbations u, w, 5P have time and space dependence as e"'*^**^^, with k real 
and n complex. We then linearize the x- and z- components of the momentum equation to obtain 

(n + C/a.)u = -— (3) 
P 

8 ()P 

{n + Ud,)w + wd,U = — - — . (4) 

P 

We expressed the perturbed continuity equation in terms of perturbed pressure, 

p pci 

where = dP/dp. We solve Equation [s] for 6P using Equation [sj and discard terms of order U/cs 
to obtain 

k2/^2/,2 (^ + U0.){-d.w) (6) 

and insert that into Equation [4] to obtain a new equation for z-momentum: 

(n + Ud,)w + wd,U = (^2^^^(ri + Ud,){d,w)^ (7) 
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which can now be written as a differential equation for w (taking U and Cs constant throughout 
the post-shock layer), 

^ >9z + ,o ""o/o gf -Ud,-n]w = 0. (8) 



We define 

J = y^^T^, (9) 

which describes the effective lateral wavenumber. As a wave approaches the acoustic case, = 
—k'^c'l, the wave becomes purely longitudinal and j tends toward zero. Equation [s] has the general 
solution 

w = Ae^' + Be-^' + Ce"""/^. (10) 

The system accordingly has three boundary conditions at its two interfaces: the shock and the rear 
surface. We note that the shock frame's acceleration Vg does not appear in the general form of the 
perturbations; it will enter into the system through the boundary conditions. 



The last term in Equation 10 is a consequence of the background flow U and is closely connected 
with structures convecting downstream with that velocity. It is instructive to consider the general 
solution for w in the frame of the rear surface. We introduce the coordinate z' = Ut — z. In 
addition, we will now write explicitly the implicit time-dependence e"*. The general solution 
is w = ^e("+^'^)*-^'^'+"* + 5e("-^'^)*+^'^' + Ce"^'/^. We see that the third term has no time- 
dependence in the frame of the rear layer. In the frame of the rear surface, these flow structures 
are generated by perturbations in the shock surface as the shock passes some point in space, and 
do not evolve further. Therefore, in the frame of the shock, this term describes flow structures 
convecting downstream through the flow with constant velocity U. We take the shock to have been 
perfectly planar at the instant, some time past, at which the shock's deceleration and rear layer 
formation began. This allows us to explicitly set C = at the rear layer. We assume however 
that the perturbation began sufficiently early in time that our treatment using Fourier modes is 
sufficient, so no further information from initial conditions will be incorporated at this time. 



3.2. Infinitely Thin Layer 



We recall that the dispersion relation for the thin shell instability in its most simple form. 



without the effects of compression, is in Vishniac and Ryu ( 1989 ) written in the form 



n +n c„k 



2„2t 2 



a 



(11) 



where a is the areal mass density of the (infinitely) thin layer, and all other variables are as we have 
defined them. Early work ( Vishniac|1983[ ) derived this expression for a shock of infinitesimal height 
but finite areal density. Such a shock, maintaining an infinitely thin layer height while continuing 
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to accrete mass from the incoming flow, would in our analysis be described as the limit of an infinite 
compression, — ?• 0. We should expect solutions we obtain for layers of finite thickness to approach 
Equation 11 in this limit. 



3.3. Free Rear Surface 



We construct the boundary condition describing a free layer at z = H hy applying 5P 



p{—Vs)5z at z = H, with dt6z = w. Using Equation 10 and our earlier expression for 6P, Equation 
[gJ the boundary condition becomes 



2 jVs)e^i^ + Bi-n^ 



jVs)e 



-jH 







(12) 



where C has been explicitly set to zero as discussed above. Equation 12 is a boundary condition 
well known to generate surface gravity waves, when j = k and when paired with a rigid boundary 
condition at z = 0. 

At the shock surface, we must perturb the shock momentum jump condition in the frame of the 
moving shock. The perturbed shock surface moving upward in Figure[T]sees a weaker incoming flow. 
In addition, by raising the shock surface in the hydrostatic pressure field, the effective post-shock 
pressure drops by an amount Vspdz. Our jump condition has now become 



PoiVs 



w 



pU^ + {P + V,p5z + 5P), 



:i3a) 



from which we obtain a boundary condition (using poVg = pU, 5z = w/n{l — ry), and our earlier 
expression for 6P in Equation^ 




n(l — rj) 



+ 2U 



w 



0. 



(13b) 



z=0 



The expression for dtSz comes from conservation of mass across the shock. With density perturba- 
tions suppressed, as discussed above, we have a balance of mass flux with pqVs entering and pU + w 
leaving the shock, with the shock moving at speed dtSz. 



r] 



U U + w- dt5z 



implying (with dt = n) 



w 



z=0 



V. - df Sz 



(1 — r])n5z 



(14a) 
(14b) 



The third boundary condition comes from oblique shock relations. Letting /3 be the angle 
of the shock surface perturbation, continuity of the tangential flow requires to first order u ~ 
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Vsf3 = {ik)Vs§z. Applying the continuity equation of Equation [5] just downstream of the shock, and 
applying Equations [6] and 14b 



Vsk' 



n(l — r]) 



w 



z=0 



{n + Ud,)5P 
pel 



which evaluates to 



n(\ — r]) 



w 



(15a) 



(15b) 



z=0 



Simultaneously applying these three conditions (equations 12, 13b, and 15b) on w, one demands 
for nonzero solutions that the determinant of the matrix of coefficients A, B, and C, shown 
collected in Equation [TgJ must be zero, 

,2 



(n2-jV,)e^;^ 

j ^ nil-n) 

1 _ Vs 

j n{l-ri) 



(- 



-n 



jVs)e-^'^ 

Vs 







- -I- r/ -I- 
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3 



2U + 



up 



Vs 
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(16) 



n{l-ri) 

From this one obtains, with some manipulation, the dispersion relation, 

= (1 - r?)n2 + fUVs + {jVs + 2njU) x 

/ (n^ + fUVs) - {njVs + n^jU) tanh jF\ 
I {n^ + pUVs) tanh jH - {njVs + n'^jU) j ' 



(17) 



We will take, as in Vishniac and Ryu (1989), the product UVg to be equivalent to an average sound 
speed squared (c^), which we shall not henceforth distinguish from the sound speed of material 



compressibility. The qualitative classification of solutions to Equation 17 depends strongly on the 
layer thickness H, specifically on its relation to the compressible scale height [714/1^41 = Cg/|T4|- 
We shall explore this dependence in what follows. 



We will now investigate the range in which wavelengths of perturbations are not much shorter 
than H, and will approximate tanhjiJ ~ jH. The existence of the critical H is easiest to see in the 
limit of very strong, highly compressive shocks (C/ — )• while — )• oo in such a way that UVg = 
and Vs remain constant). By expanding j, we may write the dispersion relation as, 

V? ' 



-SZ - k^V?S = 



(18a) 



where we have introduced scale factors 

T = 2-r] 



S 

z 



1 + 



1 + 



H 

2r?(2-77)(5-l) 
S 



(18b) 
(18c) 
(18d) 
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For strong shocks, T ~ (7 + 3)/(7 + 1), in which any effects of strong radiation are included in 7 



as an effective polytropic index describing the total density increase at the shock ( Liang and Keilty 



2000). Z is typically close to 1. Solutions of Equation 18 shown in Figure [2| yield instability for 



k in the range ki < k < k2, centered around a wavenumber of maximum instability km, where 

\V, ' 



h 

km. 



-S 





c2 


\Vs 








\Vs 







2T- Z- 2\Jt'^ - TZ 
2T- Z + 2\/r2 - TZ 
f 



(19a) 
(19b) 
(19c) 



We find that ki and /c2 are real for 5 < 0, requiring Vg < and H < c1/\Vs\, conditions defining 
a decelerating shock and a layer width shorter than a scale height. The width of the region of 
instability k2 — ki may be evaluated and is greater than zero when 



H <^ 



cl (1 - 5r/ + 2rf 



(20) 



This provides either an additional upper bound for layer width H or lower bound for compression 
1/t/, given one of the two parameters. This requirement becomes unfulfillable if the right-hand 
side of Equation 20 becomes zero, necessitating a minimum compression of approximately 4.6 to 
produce the instability. If these conditions cannot be satisfied, the solutions become stable for all 
k. 

For the high compression limit T = 2, Z = 1, the critical wavenumbers take the values 



k2 

km 



-s 



3 - \/8 - 0.293fe„ 



-s 



3 + \/8 ~ i.imkr, 



V2 



(21a) 
(21b) 
(21c) 



The solutions for growth rate at the fastest growing wavelength are 



(22) 



which allow us to verify that for k 
assumption. 



km, jH = (A;^„ + n'^/c'^)H remains small, validating our 



In the opposite limit of shock strength, as rj = I and the shock is removed from the system, 
T = 1 and the bending waves asymptotically approach unmodified acoustic waves at high k. 



Fig. 2. — Plot showing solutions of Equation 18. uj = Im(n) vs. k. The dashed line denotes the 
region of instability, where Re(n) is nonzero. 




Fig. 3. — Plot showing solutions of Equation 24. u = Im(n) vs. k. The dashed line denotes the 
region of instability, k < kcr, where Re(n) is nonzero. 
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3.4. Limiting Behavior of Solutions 

To investigate individually the effect of various terms, we may make several further simplifying 

0, as well as 

fc^. We then obtain from Equation 



assumptions to Equation 17 We will consider both the high compression limit, rj 
the limit of negligible compressibility, 



17 



the equation 



= + A;2c2 + {kVs + 2nkU) x 

(n^ + k^UVs) - {nkVs + n^Uk)tanhkH\ 
(n3 + k'^UVs) tanh jH - {nkVs + n'^Uk) I ' 



(23) 



The 2nUk term in Equation 23 which stems from the same physical source as the term 



discarded in Equation 6b of Vishniac and Ryu (1989), contributes to damping and shock stability 



in the high k limit. It was demonstrated in early work, such as that by Freeman (1955), that we 



expect stability for shocks separating two simple spaces of homogenous material. Accordingly, in 
systems with decelerating shock-bounded dense layers, as we tend to wavelengths short compared 
to the width of the layer, the dynamics must approach this stable limit ( Vishniac|1995 ). The correct 



rate of damping is however beyond the scope of our assumptions. Ishizaki and Nishimura (1997) 



have shown that the acoustic modes within the shocked material, which we have suppressed, play 
a role in stabilizing the shock. 



The limit of an indefinitely thin layer is approached, in the notation of Equation |23[ by taking 
the limit of negligible post-shock flow [/ — )■ 0, rearranging the dispersion relation as 



4 I 2/ 2 2 7 2t>2 

n + n k — k V„ 



1 + 



H 



0. 



(24) 



This shows that, in these limits, we regain the form of the Vishniac dispersion relation (Equation 
11). We also see that, for H less than the scale height and Vg < 0, the quantity in square brackets 
becomes negative, while this quantity is positive for large H or positive Vs- This means that the 



solutions to Equation 24 have the signature of the Vishniac thin layer dispersion for an accelerating 
shock, except when H lies within a scale height for a decelerating shock, H < —cl/Vs- Solutions 
when H is in that range appear as shown in Figure [3j 

The region of instability is A; < kcr, with a maximum growth at km, where 



|V:J^/=5 



k — 9- 



(25a) 
(25b) 



Compared with Equations 21 and Figure [2| we see that the principal result of removing the effects 
of compressibility is to eliminate the region of stability near /c = 0. We also see that the bending 
modes now travel asymptotically for high k with the full speed of sound, where previously they 



Fig. 4. — Numerical solutions of Equation 17, showing flow patterns of the perturbation {u,w) 
within the layer and relative phase of surface perturbations, with (a) -ff = 50 • 10~^ m, (b) H = 
110-10-6 m, (c) H = 190-10-6 m, for a shock system with Vs = 120-10^ m/sec, Vs = -5-10^^ m/sec^ 
r] = 0.05, displaying in each case a perturbation with k = 5210 m-^. 



moved at cj/VT. Lastly, all thicknesses H less than a scale height now produce instability, rather 



than following the more restrictive condition in Equation 20 



We note that for if << c^/ 



24 



becomes very large. As H 



\Vs\, the rightmost term in Equation 
becomes very close to zero, one perhaps expects this term to level off at the value in Equation 11 
we will explore this limit below. 



4. Post-Shock Flow Patterns 

Figure |4] shows a numerical solution of Equation 23 for a shock system with three different 
thicknesses. The shock system has a scale height c^/|ys| of 144- 10-^ m. One can see that for 
the very thin layer in Fig|4]^a), the flow pattern is most similar to that of a surface wave. As 
the post-shock layer increases in thickness through Figs |4]^b) and (c) , the flow pattern evolves to 
contain vorticity features. We speculate that the transition at the scale height corresponds to a 
layer thickness in which a complete cell is localized. 
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We remark that in the numerical solution of Equation 23 we find that the shock and rear 
surfaces' perturbations achieve different phase. Since the fluid inside the layer is constant in density, 
this will lead to a corresponding perturbation of areal density of the layer that might be observed. 
The physical connection is therefore maintained with the theory described by Vishniac (1983), in 
which dynamics causing variation in areal density of the post-shock layer leads to overstability in 
the shock. These plots may be compared with Figures 7-10 of Bertschinger (1986), which show 
similar vortical structure, though without boundary phase shifting. 



Further Considerations and Conclusions 



5.1. Connections to the Infinitely Thin System 

We have seen that the characteristic fourth-order nature of the Vishniac instability, as derived 



in Equation 17 follows from allowing perturbations on both surfaces of the post-shock layer. We 
note that while the Vishniac derivations contain an instability source in the product V^Pj/o", our 

■ contains a source term V?- This difference follows from Vishniac's 



17 



dispersion relation in Equation 
assumption that the post-shock layer is thin and that the difference between thermal backing 
pressure and ram pressure together with geometric factors (such as spherical divergence of the 
shock) are the fundamental sources of the deceleration. We have instead worked with planar 
shocks and assumed deceleration to stem primarily from mass accumulation and energy loss from 
the system, for example by strong radiative cooling, and a hydrostatic distribution within the layer 
to be the dominant contributor to pressure variation. 



Despite these differences in approach, we can in fact derive Equation 11 from Equation 24 
immediately. We identify the sound speed at the shock surface with local post-shock fluid variables 



PjO) ^ P^ - pVsH 
P P 



(26) 



We have implicity set the polytropic index 7 = 1, which is consistent with our assumption in 
Equation [24] that we are in the infinitely compressive limit ?? = 0. However, we do not expect 
Equation 26 to be in general consistent with our other definitions of Cg, except in the limit of an 
infinitely thin shell, — )• 0. Keeping this in mind, we see that inserting Equation 26 and a = pH 



into the term in square brackets in Equation 24 one obtains Equation 1 1 Our derivation therefore 



is found to agree with the earlier results of Vishniac in the appropriate limits. 

We comment on the different solutions to Equations [TT| [Tsj and [24j The range of stable 



solutions when H < cj/\Vs\ but is insufficiently thin to satisfy Equation 20, is unique to the 



20 



compressible, finite-shock strength solutions of Equation [T7| The oscillating instability which exists 
when Vs < and H satisfies Equation 



is the case of interest in which collective modulation of the 
boundary layers results in the growth of structure. The non-oscillating instability which appears 
when > is recognized as the Rayleigh- Taylor instability of the rear layer under acceleration. 
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The non-oscillating solutions of Equations 18 and 24 when Vg < but H > c^/iygl are of a 



different nature than the other cases. The system under perturbation was constructed by equating 
the pressure P immediately behind the shock with the ram pressure of the incoming material. The 
pressure profile then decreased hydrostatically with distance from the shock. When H exceeds a 
scale height, the most distant pressures obtained in this fashion become negative. The "instability" 
in this case is a response of the system to inconsistent initial conditions. In Vishniac's equation, 
this corresponds to the case where one assigns Vg, Pi < 0. 



Compared to Equation 11 , Equations 17 and |24| have the property of being written in terms of 
the rear layer height and variables defined locally at the compression front, with few assumptions 



regarding the structure throughout the layer, while Equation 11 is properly understood as dealing 
with quantities averaged over the layer height. This difference allows one to straightforwardly 



identify from Equation 24 the combination of system variables which lead to the transition at the 



scale height. Equation 17 features the same behavior extended to general post-shock U and finite 
rj, with appropriate corrections leading to transition at a fraction of the scale height. We expect the 
constant density solution to be applicable within a scale height, beyond which modeling the layer 
as a region of constant density will not be as appropriate as an exponential or self-similar profile. 



5.2. Experimental Observations 

We conclude with some discussion of experiments featuring strongly decelerating planar shocks. 
Experiments which intend to reproduce this instability must feature sufficient lateral space for the 
growing perturbations. Very early in the experiment's evolution, the post-shock layer thickness 
will be necessarily small, and we assume H « c^/\Vs\- We see from the results of the proceeding 



compressible analysis (Equations 21 and 22) that to allow maximum growth one must afford the 



experiment lateral dimensions A > 2ttCs\/ H/2\Vs\, where if is a characteristic or average layer 



thickness of the system. The evolution will occur within a growth time scale t = y 8H/\Vs\. 
Conversely, if one wishes to eliminate entirely this instability one should construct an experiment 



with lateral dimensions A < 2.6 Cs\/ H/\Vs\. For the experiments discussed above by 



Reighard et al. 



(2006), the values of preferred minimum distance and time are approximately 400 - 500 /um and 9 



13 ns, conditions which are achievable by the reported experiment. 
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